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ABSTRACT 

We report the detection of the lio — In ground state transition of ortho-H2D+ 
at 372.421 GHz in emission from the young stellar object NGC 1333 IRAS 4A. 
Detailed excitation models with a power-law temperature and density structure 
yield a beam-averaged H2D+ abundance of 3 x 10^^^ with an uncertainty of 
a factor of two. The line was not detected toward W 33A, GL 2591, and 
NGC 2264 IRS, in the latter source at a level which is 3 — 8 times lower than 
previous observations. The H2D+ data provide direct evidence in support of 
low-temperature chemical models in which H2D+ is enhanced by the reaction 
of H^ and HD. The H2D+ enhancement toward NGC 1333 IRAS 4A is also 
reflected in the high DCO^/HCO^ abundance ratio. Simultaneous observations 
of the N2H+ 4 — 3 line show that its abundance is about 50 — 100 times 
lower in NGC 1333 IRAS 4A than in the other sources, suggesting significant 
depletion of N2. The N2H+ data provide independent lower limits on the H.^ 
abundance which are consistent with the abundances derived from II2D+. The 
corresponding limits on the H^ column density agree with recent near-infrared 
absorption measurements of HjJ' toward W 33A and GL 2591. 

Subject headings: ISM: abundances — ISM: molecules — molecular processes — 
radio lines: ISM — submillimeter 
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1. Introduction 

The recent detection of tlie Hg" ion in interstellar clouds through its infrared vibration- 
rotation lines (Geballe & Oka 1996; McCall et al. 1998) is an important confirmation of 
the gas-phase chemical networks (Herbst & Klemperer 1973; Watson 1973). Because of its 
symmetry, Hj^ has no allowed rotational transitions contrary to its deuterated isotopomer 
H2D"'" which has a large permanent dipole moment (Dalgarno et al. 1973). Thus, H2D"'" 
is important as a tracer of B.^ with transitions that can be searched for in emission. 
In addition, it is widely believed to play a pivotal role in the interstellar ion-molecule 
chemistry at low temperatures where significant enhancement of deuterated molecules 
occurs as a result of fractionation (e.g. Herbst 1982; Millar et al. 1989). This process is 
initiated by the isotope exchange equilibrium reaction Hj + HD ^ H2D+ + H2 (1) which 
is shifted in the forward direction at low temperatures (Smith et al. 1982; Herbst 1982). 
The formation of H2D"'' is followed by deuterium transfer reactions with e.g. CO to form 
DCO^, and the H2D+ enhancement is reflected in the observed large abundance ratios of 
e.g. DCO+/HCO+, NH2D/NH3, and DCN/HCN in cold clouds (e.g. Wootten 1987; Butner 
et al. 1995; Williams et al. 1998). 

Over the last 20 years numerous attempts have been made to detect the 1 10 — In ortho- 
H2D"'" and Iqi — Oqo para-H2D+ ground state lines at 372 and 1370 GHz, respectively. These 
searches have mainly been done with the Kuiper Airborne Observatory (KAO) (Phillips et 
al. 1985; Pagani et al. 1992a; Boreiko & Betz 1993), and a possible absorption feature at 
1370 GHz has been reported by Boreiko & Betz (1993) toward Orion. Observations from 
the ground are very difficult since the 372 GHz line is at the edge of a strong atmospheric 
water absorption line, while the atmosphere at 1370 GHz is almost completely opaque. 
With the advent of new submillimeter receivers equipped with sensitive niobium SIS mixers 
it has become possible to search for weak ortho-H2D"'" lines from high, dry sites such as 
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Mauna Kea. Indeed, a ground-based search for this hne from the Caltech Submillimeter 
Observatory (CSO) by van Dishoeck et al. (1992) yielded hmits which are up to a factor 
of hundred more sensitive than those obtained with the KAO. Comparison with chemical 
models suggested that only a factor of a few improvement would be needed to detect the 
line. With the new facility receiver RxB3 at the James Clerk Maxwell Telescope (JCMT)0 
such an improvement in sensitivity is now achievable. Here we report the detection of the 
H2D+ 372.421 GHz line toward NGC 1333 IRAS 4A and significant upper limits toward 
W 33A, GL 2591, and NGC 2264 IRS. Simultaneous observations of the N2H+ 4-3 line 
at 372.672 GHz toward these young stellar objects (YSOs) are used to place additional 
constraints on the H;^ abundance. 

2. Observations 

The observations of the lio — In ground state transition of ortho-H2D+ at 372.42134 
GHz (Bogey et al. 1984) were done with the JCMT on August 31, September 15 and 18, 
1998 during three nights of very good submillimeter transparency with a zenith opacity 
at 225 GHz below 0.05. The dual-polarization heterodyne receiver RxB3 was used (Avery 
1998). Both mixers were tuned to 372.5469 GHz in the upper sideband. The big advantage 
of RxB3 is that it has a dual-beam interferometer which allows single-sideband (SSB) 
operation, enhancing the sensitivity and calibration at 372 GHz considerably. The digital 
autocorrelator spectrometer (DAS) was split into four parts of 125 MHz. This setup allows 
observations of both lines in two orthogonal polarizations simultaneously with a spectral 



^ The JCMT is operated by the Joint Astronomy Centre in Hilo, Hawaii on behalf of the 
Particle Physics and Astronomy Research Council in the United Kingdom, the Netherlands 
Organisation for Scientific Research, and the National Research Council of Canada. 
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resolution of 376 kHz (= 0.3 km s^^ at 372 GHz). Typical SSB system temperatures 
including atmospheric losses were about 1200 K. The effective total integration time was 
7 hours on NGC 1333 observed over two nights, 2.7 hours on W 33A, 4.3 hours on GL 
2591, and 4 hours on NGC 2264. The absolute calibration uncertainty is estimated at 30%, 
and the relative calibration between the H2D+ and N2H+ lines is much better. The JCMT 
beamsize at 372 GHz is 13" FWHM, the main beam efficiency is 57%. JCMT data on 
H^^CO^ and DCO+ were taken from the literature (see below), except for DCO+ toward 
W 33A and GL 2591, for which the 3 — 2 transition at 216.113 GHz was observed with 
receiver RxA3. The beamsize at this frequency is 21" FWHM and the main beam efficiency 
is 70%. 

The observed H2D+ and N2H+ spectra are presented in Figure 1. The source and line 
parameters are listed in Table 1. The H2D+ line is clearly detected with T4 = 0.08 ± 0.03 
K toward NGC 1333 IRAS 4A, and is seen in spectra of both nights. The velocity width 
shows good agreement with the N2H"'" line width while the velocities are offset by about 0.5 
km s^^. Comparison with the line survey of Blake et al. (1995) shows that such an offset 
is small and common for this region. No H2D+ emission was detected toward NGC 2264 
IRS, W 33A, and GL 2591. Assuming the same width as the N2H+ line, 2a upper limits of 
Tl < 0.02 - 0.04 K are obtained. For NGC 2264 this limit is about a factor of eight below 
the possible feature of Phillips et al. (1985), and a factor of three below the limit reached 
by van Dishoeck et al. (1992). Note that the N2H"^ emission toward NGC 1333 is much 
weaker than that toward the other sources. No other lines were detected in the 125 MHz 
bands. 
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3. Analysis 

Model calculations were performed to determine the abundances of H2D"'", N2H"'", 
HCO"*", and DCO^ using a power-law density structure n = no{r/Ro)~°', as described in 
van der Tak et al. (1999). In these models the radial dust temperature profile is calculated 
from the observed luminosity and no is determined from submillimeter photometry, which 
probes the total dust mass. The grain heating and cooling are solved self-consistently as a 
function of radius, r, using grain properties from Ossenkopf & Henning (1994). The outer 
radius (Ro) is determined from high resolution submillimeter line and continuum maps. 
The exponent a is constrained by modeling the relative strength of emission lines of CS 
and H2CO at the central position over a large range of critical densities with a Monte Carlo 
radiative transfer program, assuming Tk = T^nst- Data were taken from Blake et al. (1995) 
(NGC 1333 IRAS 4A), de Boisanger et al. (1996) and Schreyer et al. (1997) (NGC 2264 
IRS), and van der Tak et al. (1999 and in preparation) (GL 2591 and W 33A). For NGC 
1333 IRAS 4A, where CS is heavily depleted, a = 2 was taken based on the analysis of the 
continuum visibilities in interferometer data by Looney (1998). 

Given the calculated temperature and density structure, the radiative transfer models 
were run to determine the abundances, assuming initially a constant abundance throughout 
the envelope. Both the ortho-Il2D+ and para-H2D"'" ladders have been considered since 
their spin states are coupled through reactive collisions with H2; thus the para Oqo level is 
the true rotational ground state. A de-excitation rate coefficient of 1.0 x 10^^*^ cm^ s^^ has 
been used for all inter-ladder transitions. See Herbst (1982) and Pagani et al. (1992b) for a 
detailed study of the ortho/para ratio. The lower level of the lio — In transition lies at 86 
K. The excitation energy of the lio level is 18 K relative to the In level, and the critical 
density for this transition is about 2 x 10^ cm~^. 

The calculated abundances are listed in Table 2. Toward NGC 1333 we infer a 
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beam-averaged abundance a:(H2D+) = 3 x 10^^^. Upper limits on the abundance toward 
NGC 2264, W 33A, and GL 2591 are < 1 x lO^^^ The N2H+ abundance ranges between 
10^^^ toward NGC 1333 and 10^^ toward W 33A. All derived abundances have an absolute 
uncertainty of a factor of two, due to the uncertainties in the dust opacities and CO 
abundances. The relatively high N2H+ abundance toward NGC 2264 was already found 
by van Dishoeck et al. (1992), who noted that nearly all of the gas phase nitrogen must 
be in the form of N2 in this cloud. Since N2H+ is formed mainly by the reaction of H3" 
and N2, the latter observations provide an independent lower limit on the H3" abundance. 
Destruction occurs mainly via reactions with CO, O, and electrons. Considering CO 
destruction only, ^(H^) ^ 0.5 n(N2H+) a;(C0)/a;(N2) (2). Assuming 50% of the nitrogen is 
in N2, a;(N2) = 5 x 10'^6(N2), a;(CO) = 2 x 10-^(5(CO) and equal amounts of depletion S 
for CO and N2, this yields ^(Hjj") ^ 2a;(N2H+). These limits are listed in Table 2, and are 
consistent with the upper limits derived from the H2D+ observations using a theoretical 
H2D+/H^ ratio, see §4. 

4. H2D+/H;^ chemistry 

The above analysis assumes constant abundances throughout the YSO envelopes. 
In reality, the H2D"'" abundance is a strong function of temperature and position. In 
chemical equilibrium, the H3" abundance can be written as ^(Ilg") = (/Tikxn(X) with 
n(X) = n(H2)a;(X). X refers to any of O, C, CO, O2, N2, H2O, ... which are the principal 
removal agents of B.^ via the proton transfer reactions Hf + X ^ XH^ + H2 (3), where kx 
are the rate coefficients (taken from the UMIST database, see e.g. Millar et al. 1997) and ( 
is the cosmic-ray ionization rate (taken to be 5 x 10^^^ s^^). A simple chemical model for 
the formation and destruction of 1120^^ yields 

a;(H2D+) x{BI))kf + x{D)kD 



a;(H^) x{e)ke + Ekxx{X)+kr' 



(4) 



where kf and kr are the forward and backward rate coefficients of reaction (1), kn is the rate 
coefficient for formation of H2D+ through the reaction H^ + D, and ke is the rate coefficient 
of the electron recombination of H2D+ (see e.g. Caselh et al. 1998 for a compilation of 
values). We assumed x{WD) = 10a;(D) = 2.8 x 10~^ throughout. 

The above H;^ and H2D+ chemical equations were included in the power-law models, 
and abundances at each position were calculated for the appropriate temperature and 
density. We have fixed the expression for /c^ at T < 20 K to its value at 20 K, to ensure 
that x(H2D+) < a;(Hj) throughout. For simplicity, only X=CO was considered and the 
electron recombination was neglected. The CO depletion is inferred from C^^O observations 
as described in van der Tak et al. (1999 and in preparation), and are listed in Table 2. For 
a homogeneous temperature and density structure, our model agrees well with the models 
of Millar et al. (1989) and Pagani et al. (1992b). 

The power-law model results for NGC 1333 IRAS 4A are presented in Figure 2. Using 
these abundances, the H2D'^ emission has been calculated, most of which originates from 
gas at T = 25 — 35 K. The model intensity agrees within 30% with that measured toward 
NGC 1333 IRAS 4A, and is consistent with the upper limit toward GL2591. They are a 
factor of 2 — 4 larger than the upper limits toward NGC 2264 IRS and W 33A, respectively. 
Most likely, this small discrepancy results from effective removal of H;^ by other species 
than CO and/or by an overestimate of the size of the cloud. The II3" column densities 
computed for W 33 A and GL 2591 agree within a factor of 2 — 3 with the directly observed 
column densities by Geballe & Oka (1996) and McCall et al. (1998). This good agreement 
between models and observations provides the strongest support for the basic H3" and H2D+ 
chemical networks in dense, cold clouds. 

The simple chemical networks described above (Eq. (2) and (3)) have also been used 
to model the N2H+ abundance from the N2H+/H2D+ ratio. Assuming equal amounts of 
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depletion for CO and N2, agreement within a factor of two between the measured and 
modeled N2H"'" hne intensities is obtained. The best fitting N2 abundances are included in 
Table 2. 

The derived H2D+/HJ]' ratios are also compared with the observed DCO+/HCO^ 
ratios in Table 2. The enhancement in H2D+/H^ is clearly reflected in the DCO^/HCO+ 
ratio, as expected since the latter species are directly formed by reactions of the former 
with CO. For NGC 1333 IRAS 4A, DCO+/HCO+ = 0.5 H2D+/HI, where the factor of 
0.5 is consistent with the statistical branching ratio of 1/3 within the uncertainties. The 
DCO+/HCO+ ratio toward NGC 1333 IRAS 4A is a factor of 5 - 50 higher than that toward 
the other sources. This can be explained by the difference in physical structure. In cold, 
dense (pre-)protostellar cores like NGC 1333 IRAS 4A where the CO and N2 depletions are 
extreme, the II2D+ abundance is enhanced, because of the low temperature and because 
the main removal reactions of II3" are suppressed. The H2D'*' abundance will increase even 
stronger than that of H^ since reaction (1) becomes the main destruction channel of H^. 
In the case of NGC 2264, W 33A, and GL 2591, where the temperatures are higher and the 
depletion of CO and N2 is less, reaction (3) becomes the main removal path of Hg". 

The difference in physical structure may have its origin in the different stages of 
protostellar evolution. In particular, NGC 1333 IRAS 4A has been classified as a Class 
object (Andre & Montmerle 1994), and is thus in a very early stage, with a large spatial 
separation between the quiescent and shocked regions (Blake et al. 1995). 

5. Conclusions 

With the current sensitivity of heterodyne receivers it is now possible to study the 
ortho-H2D"'' lio — In 372.421 GHz line profile in emission in the early stages of star 
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formation deep inside dense molecular clouds. Its importance lies in the fact that it is a 
tracer of H^ and that it provides information on the deuterium abundance and temperature 
history of a cloud and on the chemical evolution during star formation. Further observations 
of the H2D+ and H^ lines in a sample of very young Class and Class I YSOs will therefore 
be very valuable. 

Since the ortho-H2D"'" ground state is at 86 K, the lio — In line traces both the 
warm and cold regions, although the H2D"'" enhancement will be strongest in the coldest 
regions. Observations of the Iqi — Oqo para-H2D+ ground state line at 1.37 THz toward 
the continuum of embedded YSOs may reveal cold H2D''" in absorption. The dual 
channel German REceiver for Astronomy at THz frequencies (GREAT) to be flown on the 
Stratospheric Observatory For Infrared Astronomy (SOFIA) would allow such observations. 
Combined with 372 GHz observations from the ground, the total abundance and the 
relative population of the ortho- and para-modifications may be determined which provides 
information on the formation, destruction, and excitation processes. Simultaneous deep 
observations of the HD J = 1 — (2.7 THz) and para-H2D'^ ground state lines toward 
YSOs may yield a direct measure of the (variation in) H^ abundance over the cloud, and 
thus of the cosmic-ray ionization rate. 

It is a pleasure to thank Lome Avery for discussions and test observations at 372 GHz, 
Henry Matthews for the non-standard setup for the backend, and Leslie Looney and Lee 
Mundy for sharing their model results on NGC 1333 prior to publication. The observations 
would not have been possible without the fiexible observing strategy of Remo Tilanus. 
Gerd-Jan van Zadelhoff, Fred Baas, and Jane Greaves did an excellent job in carrying out 
these observations in service. This work is supported by NWO grant 614.41.003. 
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Fig. 1.— Observed spectra of H2D+ lio - In at 372.421 GHz and N2H+ 4 - 3 at 372.672 
GHz. The N2H"'' spectra have a resolution of 0.3 km s~^. The H2D+ spectra have been 
smoothed to 0.6 km s~^. 

Fig. 2. — Power-law model abundances of H^, H2D+, and N2H+ as function of density and 
temperature. 



-14- 



Table 1. Observations'^ 



Source'' Molecule Transition TJ AV Vlsr 



(K) (km s-i) (km s-l) 



N1333 


H2D+ 


lio — 111 


0.08(0.03) 


1.2 ±0.3 


7.4 ±0.2 




N2H+ 


4-3 


2.57(0.03) 


1.35 


6.94 


N2264 


H2D+ 


lio - 111 


< 0.02== 








N2H+ 


4-3 


4.51(0.03) 


2.67 


8.01 


W33A 


H2D+ 


lio — 111 


< 0.04== 








N2H+ 


4-3 


3.06(0.06) 


4.62 


37.40 


GL2591 


H2D+ 


lio - 111 


< 0.02== 








N2H+ 


4-3 


1.41(0.04) 


2.89 


-5.82 



''Values in parentheses represent la statistical uncertainties. The absolute 
uncertainty of the intensity is 30%; AV and Vlsr are accurate to better than 
0.1 km s^i. 

''Position (B1950): NGC 1333 IRAS4A: a = 03''26™04!8, <5 = 
+31°03'14"; NGC 2264 IRS: a = 06''38™25?0, <5 = +09°32'29"; W33A: a = 
18im™44!2, <5 = -17°52'56"; GL 2591: a = 20''27™35!8; <5 = +40°01'14" 

'=20- limits: NGC 2264 IRS: AV = 2.5, GL 2591: AV = 3.0, and W33A: 
AV = 4.5 kms-i. 



-15- 



Table 2. Excitation model parameters and deduced abundances'* 



Source a no T 

Ri, Ro 



(cm-3) (K) H2D^ 



Molecular Abundances 



N2H+ 



jj+b 



TT+C HaP^ 

"3 „+ 



DCO+ 
HCO+'^ 



co<= 



N2 



N1333 2 1.7(6) 

N2264 1.5 1.5(4) 

W33A 1 2.1(4) 

GL2591 1.25 3.5(4) 



318,13 3(-12) 

293,18 < l(-ll) 

280,16 < l(-ll) 

350,30 < l(-ll) 



l(-ll) 


2(-10) 


>2(- 


-11) 


2(-2) 


l(-2) 


4(-6) 


4(-7) 


6(-10) 


< 3(- 9) 


>1(- 


- 9) 


< 3(-3) 


3(-3) 


6(-5) 


4(-5) 


1(- 9) 


< 4(- 9) 


>2(- 


- 9) 


< 3(-3) 


l(-4) 


5(-5) 


2(-5) 


5(-10) 


< 3(- 9) 


>1(- 


- 9) 


< 3(-3) 


4(-4) 


2(-4) 


9(-5) 



^From statistical equilibrium calculations using the appropriate temperature and density structure as function of distance r 
to the YSO, n{r) = Uoir/Ro)^" where Ro is the outer radius of the model envelope: NGC 1333 IRAS 4A: 3.1(3) AU, NGC 
2264 IRS: 4.7(4) AU, W33A: 2.4(5) AU, GL 2591: 3.1(4) AU, and R, = Ro/300 is the inner radius. Notation a(b) indicates 
a X 10*". The accuracy of the deduced abundances is a factor of two. 

''From H2D+ using a theoretical Hg /H2D+ ratio at the effective temperature from which most of the emission arises (see 
Figure 2). 

'^From N2H+ analysis (see text). 

'iFrom Hi3cO+ assuming HCO+/Hi3cO+ = 60. 

''From C^'^O assuming CO/C'^O=2500 and using the appropriate A'^(H2) from submillimeter dust emission in a 13" beam: 
NGC 1333 IRAS 4A: 3.1(23) cm-2, NGC 2264 IRS: 1.2(23) cm-2, w 33A: 5.2(23) cm-2, QL 2591: 1.3(23) cm-2. 
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